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ABSTRACT 

Approximation  of  the  differential  equation  of  the  single 
degree  of  freedom,  damped  oscillator  by  the  first  order  central  finite 
difference  equations  and  two  Hermitian  finite  difference  equations  is 
investigated.  An  error  analysis  is  made  between  the  solution  to  the 
differential  equation  and  the  solutions  to  the  finite  difference 
approximations  for  various  values  of  damping.  The  results  of  the  error 
analysis  indicate  that  the  Hermitian  approximations  have  less  error 
than  the  first  order  difference  approximation  for  the  same  size  of 
increment.  Furthermore,  the  employment  of  the  Hermitian  method  does 
not  materially  increase  the  execution  time  on  the  digital  computer  over 
that  required  by  the  first  order  difference  equations.  Thus  the 
Hermitian  equations  are  superior  to  the  first  order  finite  difference 
equations  for  the  damped  oscillator  problem. 
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CHAPTER  I 
INTRODUCTION 

In  the  quest  for  accurate  numerical  solutions  to  differential 
equations,  considerable  effort  has  been  devoted  to  the  determination 
of  the  errors  associated  with  various  finite  difference  schemes.   This 
work  has  gained  importance  with  the  recent  widespread  use  of  the 
digital  computer.  Due  to  the  high  value  placed  on  computer  time  and 
to  the  limitations  on  core  storage  space  within  the  computer,  the 
accuracy  and  efficiency  of  the  numerical  method  employed  are  of  the 
utmost  importance. 

In  general,  many  engineering  problems  can  be  classified  as 

2 
either  boundary  value  problems  or  initial  value  problems.   For  the 

numerical  solution  to  boundary  value  problems,  the  most  common  method 

used  is  a  form  of  the  central  finite  difference  operators.   Several 

forms  of  these  operators  can  be  found  in  Table  III  of  the  Appendix  to 

reference  1.   Some  of  these  are  frequently  used  to  solve  engineering 

problems  while  others  seem  either  to  be  impractical  or  to  appear 

seldom  in  current  literature.   The  Hermitian  equations,  which  are 

included  in  this  table,  fall  into  the  last  category.   For  initial 

value  problems,  numerical  methods  of  solution  such  as  Adams,  Runge- 

Kutta,  and  other  forms  of  the  predictor-corrector  methods  are  usually 

employed.  The  Hermitian  method  is  also  applicable  to  initial  value 

problems  but  is  seldom  used. 

The  purpose  of  this  paper  is  to  evaluate  the  merits  of  the 

Hermitian  method  when  applied  to  a  second  order,  ordinary,  linear 

differential  equation.   This  evaluation  will  be  made  by  comparing  the 

solution  to  the  exact  differential  equation  with  the  solution  to  two 

Hermitian  finite  difference  approximations  and  with  the  solution  to 

the  first  order  central  finite  difference  equation.   In  this  manner, 

the  difference  in  accuracy  between  these  finite  difference  methods  can 

be  determined,  and  the  merits  of  the  Hermitian  method  established. 
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The  equation  chosen  for  this  analysis  is  that  of  the  single 
degree  of  freedom,  damped  oscillator.  This  equation  was  chosen 

because  it  contains  both  the  first  and  second  derivatives  and  because 

4  5 
its  exact  solution  is  well  known,  '   This  appears  to  be  the  first 


application  of  the  Hermitian  equations  to  differential  equations  con- 
taining more  than  one  degree  of  derivatives.  This  problem  is  generally 
considered  an  initial  value  problem.  However,  the  results  of  this 
analysis  are  also  applicable  to  boundary  value  problems  since  a  boundary 
value  can  be  solved  using  methods  of  solution  normally  employed  in 
initial  value  problems.  Also,  the  oscillator  may  become  a  boundary 
value  problem  if  the  displacement  or  velocity  is  specified  at  a  time 
other  than  zero.   Thus,  the  choice  of  this  example  allows  the  conclu- 
sions to  be  applied  to  both  initial  value  and  boundary  value  problems. 

The  author  wishes  to  acknowledge  a  special  thanks  to 
Professor  Robert  E.  Ball  of  the  Aeronautics  Department  of  the  Naval 
Postgraduate  School.  The  guidance  and  consultation  of  Dr.  Ball,  as 
thesis  advisor,  were  instrumental  in  the  development  of  this  paper. 
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CHAPTER  II 

DISCUSSION  OF  THE  CENTRAL  FINITE 
DIFFERENCE  METHOD  AND  ITS  ASSOCIATED  ERRORS 

There  are  three  types  of  error  associated  with  the  finite 

3 
difference  solution  to  engineering  problems.   One  type  is  called 

truncation  or  discretization  error.  The  conversion  of  a  differential 

equation  to  an  algebraic  equation  is,  in  essence,  the  approximation 

of  a  continuous  system  by  a  discrete  mode.  The  error  involved  in 

this  conversion  is  called  the  truncation  error.  A  graphic  example  of 

how  this  error  is  introduced  can  be  seen  by  examining  the  derivation 

of  the  first  order  central  finite  difference  equations  for  the  first 

and  second  derivatives  presented  in  Appendix  A. 

In  the  derivation  of  the  finite  difference  expression  for  the 
first  derivative,  only  the  first  three  terms  of  the  Taylor  series  are 
used.  The  first  four  terms  are  used  in  the  expression  for  the  second 
derivative.  The  truncation  error  can  be  approximated  by  the  first 
unused  term  in  the  Taylor  series.  Therefore,  a  decrease  in  truncation 
error  can  be  accomplished  by  using  more  terms  of  the  Taylor  series  in 
the  finite  difference  equation. 

Another  type  of  error  is  called  round-off  error.  This  error 
is  generated  within  the  computer  and  is  a  function  of  the  amount  of 
calculation  required  and  the  condition  of  the  numbers  involved.   Each 
arithmetic  computation,  such  as  addition,  multiplication,  and  expo- 
nentiation, has  a  round-off  error  associated  with  it.   The  computer 
is  limited  to  carrying  a  specified  number  of  digits.  Any  amount  over 
this  number  must  be  dropped  for  each  calculation.  To  decrease  round- 
off error,  it  is  essential  to  minimize  the  number  of  arithmetic  calcu- 
lations within  the  program  and  to  avoid  subtraction  of  large  numbers,  etc, 

The  third  type  of  error  is  called  inherent  error.  This  error, 

unlike  truncation  error  and  round -off  error,  is  not  caused  by  the  use 

3 
of  numerical  methods.   Inherent  error  is  introduced  into  the  program 

within  the  data.   If  the  data  are  not  exact,  which  is  usually  the  case, 

the  solution  cannot  be  expected  to  be  exact. 
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The  object  of  any  numerical  analysis  is  to  obtain  as  accurate 
an  answer  as  possible  without  significantly  increasing  the  computation 
time.  There  are  a  number  of  higher  order  finite  difference  methods 
available  that  have  less  truncation  error.  Hence,  they  provide  more 
accurate  solutions  than  the  first  order  central  finite  difference 
method.  Among  these  is  the  Hermitian  method.  These  methods  vary  by 
the  number  of  stations  used,  the  location  of  the  stations,  etc. 

When  central  finite  difference  methods  are  applied  to  an 
engineering  problem,  the  problem  is  reduced  to  that  of  solving  a 
matrix.  This  matrix  will  usually  be  a  band  matrix  centered  about  the 
principal  diagonal.   The  width  of  the  band  is  dependent  upon  the  number 
of  stations  included  in  the  finite  difference  equation.  Using  the  first 
order  difference  equations  listed  in  Appendix  A,  the  band  width  would 
be  three,  and  the  matrix  is  called  a  tridiagonal  matrix.   If  higher 
order  difference  equations  that  require  more  stations  are  used,  the 
band  width  widens.   By  widening  the  band  of  the  matrix,  more  computa- 
tions may  be  required  than  for  the  tridiagonal  matrix.  This  increases 
the  round-off  error  and  the  computation  time.  Furthermore,  the  tri- 
diagonal matrix  will  require  less  core  storage  than  the  wider  band 
matrix  for  the  same  element  size. 

Another  approach  used  to  decrease  the  magnitude  of  the  trun- 
cation error  involves  decreasing  the  interval  between  the  stations. 
This  is  equivalent  to  increasing  the  degree  of  the  matrix  and,  there- 
fore, introduces  greater  round-off  error.  Other  effects  of  this  pro- 
cedure are  the  increase  in  computer  storage  space  required  and  the 
increase  in  computer  running  time. 

The  Hermitian  method  considered  in  this  thesis  results  in  a 
tridiagonal  matrix  but  also  decreases  the  truncation  error.  Conse- 
quently, the  gain  in  accuracy  due  to  the  decrease  in  truncation  error 
is  not  offset  by  an  increase  in  round-off  error.  Furthermore ,  this 
method  can  be  used  to  advantage  on  those  problems  where  core  storage 
requirements  are  important. 
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CHAPTER  III 

APPLICATION  OF  THE  HERMITIAN  METHOD 
TO  THE  DIFFERENTIAL  EQUATION 


The  first  order  central  finite  difference  equations  for  the 
first  and  second  derivative,  as  derived  in  Appendix  A,  are 

yj  =2A\yJ+i  '    yJ~i   )  (3aa) 

yj  =-gr\yji-i  -2%   +    ty-,]  o.ib) 

The  Hermitian  equations  for  the  first  and  second  derivatives,  as 
derived  in  Appendix  B,  are 


(3.2a) 
(3.2b) 


where  (')  denotes  the  derivative  with  respect  to  time. 

As  shown  in  Appendix  B,  the  Hermitian  equations  have  less  trun- 
cation error  than  the  first  order  equations.  Note,  that  in  order  to 
get  this  reduction  in  trucation  error  it  was  not  necessary  to  use 
more  stations. 

Consider  the  single  degree  of  freedom,  damped  oscillator  shown 

4  5 
in  Figure  1„  The  differential  equation  of  the  oscillator  is   ' 

m  y"t  cy1  4  k  y     —   q  (3.3) 

where  m  is  the  mass  of  the  system,  c  is  the  viscous  damping  coeffi- 
cient, k  is  the  spring  constant,  and  y  is  the  displacement.   Sub- 
stituting the  conventional  first  order  central  difference  equations, 
Eqs.  3.1a  and  3.1b,  for  the  first  and  second  derivatives  into  Eq.  3.3 
leads  to 

Efty*w*fr,)+f£t\yj*ryh,)+kyj-<>    o.*> 

where  ^  has  been  replaced  with^  f   to  denote  the  time  interval. 
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FIGURE    I      THE   SINGLE    DEGREE    OF   FREEDOM, 
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Since  the  Hermitian  equations  are  not  of  the  same  form  as 
the  conventional  equations s  a  different  method  of  approximation  to 
Eq.  3.3  is  necessary.   Examination  of  the  right  side  of  Eqs.  3.2a 
and  3.2b  reveals  that  Eq.  3.2a  contains  the  first  derivative  at  the 
three  stations  j+1,  j,  and  j-1.  Likewise,  the  right  side  of  Eq.  3.2b 
contains  the  second  derivative  at  the  same  three  stations.   By  adding 
a  third  equation  expressing  y  at  j+1,  j,  and  j=l  and  modifying  one  of 
the  Hermitian  equations,  it  is  possible  to  solve  the  differential 
equations  at  these  three  stations  simultaneously.   Equations  3.5, 
3.6,  and  3.7  show  this  procedure. 


*!%..-*% *».,  )-7rk;>*,'%:, 


(3.5) 


C 

4  At 


^  -  yJ- 


/ 


,2  t>j*i  +  4y>+y 


V+& 


Adding  Eqs.    3.5,    3.6, 


gives 


(3.6) 


&W/0Vyj-/  HtfIv0^^-,    I 


(3.7) 


^  (V  2V^-< 


4A 


-h+,-yj-,  ) 


7T\w°yj*yj-, 


10  m     if       4  c        f      10  k 
— U    -f    — —  u  +  '-^-^  u, 

12     V         12    *J      12     V 


(3.8) 
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6c 
Note  that  the  right  side  of  Eq.  3„8  will  equal  zero  if  yr-  y°  .  is 

added  to  both  sides .  The  conventional  approximation  for  y" „  gives 


(3.9) 


If  Eq.  3.9  is  added  to  Eq.  3.8,  the  result  is 

m 


A 


(3.10) 


This  is  one  form  of  the  Hermitian  approximation  to  Eq.  3.3. 

To  get  the  second  form  of  the  approximation,  it  is  necessary 
to  write  Eq.  3.7  as 


6 

(3.11) 
When  Eqs.  3.5,  3.6,  and  3.11  are  added,  then 

2m 


At* 


\yjH       Hj    +yJ-i)    +2A-M#+/         Hi-*) 

(3.12) 
Note  that  the  right  side  of  Eq.  3.12  will  equal  zero  if  my'' .  is 
subtracted  from  both  sides.  From  the  conventional  approximation  for 


v1  ' 

y  y 


a?\Hh,-2Hi +  yj-t  \=m*j    (3a3) 
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If  Eq.  3.13  is  subtracted  from  Eq.  3.12,  the  result  is 

This  is  the  second  form  of  the  Hermitian  approximation  to  Eq.  3.3. 
In  this  thesis,  Eq.  3.10  will  be  called  MOD  2S  and  Eq.  3.14  will  be 
called  MOD  1. 

By  following  the  derivation  of  these  two  Hermitian  approxi- 
mations, it  should  be  clear  that  a  first  order  central  finite 
difference  equation  for  either  y" „  or  y°°.,  was  used  in  each  deriva- 
tion. For  example,  the  first  order  central  finite  difference  equa- 
tion for  y'.,  Eq.  3.9,  was  introduced  in  the  derivation  of  Eq.  3.10, 
the  MOD  2  approximation.  This  is  equivalent  to  modifying  Eq.  3.6 


to  be 
C 


%  "Hhf  )--fr[^l+l°yj+yj-l 


2Af\v+/       V-'    /       Id   VU+i  JJ     "J- 1  I      (3.15) 

Appendix  C  shows  that  the  truncation  error  in  Eq.  3.15  is  the  same 

as  the  truncation  error  in  the  first  order  central  finite  difference 

equation  for  the  first  derivative,  Eq.  3.9.  Thus,  when  the  Hermitian 

equations  are  applied  to  the  differential  equation  containing  more 

than  one  degree  of  derivative,  accuracy  can  be  gained  for  only  one 

specified  derivative  at  a  time.  The  MOD  1  approximation  has  the 

accuracy  gain  in  the  first  derivatives,  and  the  MOD  2  approximation 

has  the  accuracy  gain  in  the  second  derivative. 

Examination  of  Eqs.  3.4,  3.10,  and  3.14  reveals  that  the 

only  difference  in  these  three  equations  is  in  the  third  term. 

Equations  3.10  and  3.14  allow  values  for  y.,1  and  y„  .,  while  Eq.  3.4 

j+l      j-i 

does  not.  This  ability  to  let  y  vary  over  the  interval  j+l  to  j-1  is 
the  reason  truncation  error  is  decreased.  If  the  variation  of  y 
happens  to  be  linear  over  the  three  successive  stations  j+l,  j,  and 
j-1,  Eqs.  3.10  and  3.14  then  both  reduce  to  the  conventional  equation, 
Eq.  3.4. 
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CHAPTER  IV 
COMPARISON  OF  THE  SOLUTIONS 

The  three  finite  difference  approximations  to  the  differen- 
tial equation,  Eq.  3,3,  are  repeated  here  for  convenience. 
Conventional 


m 


At' 


*J-> 


(3.4) 


MOD   2 


m 


At; 


yJ«~2Hj  ^-i  1^   23fl^'^-'   ) 


MOD    1 


/2  It/-*-/        #     ty-/  I 


(3.10) 


The  exact  solution  to  the  differential  equation  is  of  the 


form 


y  =  A   e 


at 


(4.1) 


where  A  is  an  arbitrary  constant.   Substituting  Eq.  4.1  into  the 
differential  equation,  Eq.  3.3,  leads  to 


met +  ca      +■    k 


Acaf        =  0 


If  Eq.  4.1  is  to  be  valid  for  all  values  of  t,  then 

m  a ■+■  c  a  +  k    =    0 


Thus, 


a  = 


[-c  ±   yc2  -4mk 
2  m 


(4.2) 
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If  CO  ~—\     I     K  and      Y   —  •' —        ,    then  Eq.    4.2  becomes 

V  m  b    2mu) 

a    =   —  to  i     ±    I  to  -\j  i  ^r  2 

b  (4.3) 

K     is  the  nondimensional  damping  coefficient.  A  value  of 
L   less  than  1.0  means  the  system  is  underdamped,  and  the  motion  is 
oscillatory.  A  value  of  K    greater  than  1.0  means  the  system  is 
overdampeds  and  the  motion  is  non-oscillatory.  When  Y      is  exactly 
equal  to  1.0,  the  system  is  critically  damped.  For  this  analysis, 

the  special  case  of  critical  damping  is  not  considered. 

2 

The  solution  assumed  for  the  finite  difference  equations  is 

y,  =  a  \J 

(4.4) 
where  j  is  the  station  number.   Solutions  for  )^  that  satisfy 
Eqs.  3.4,  3.10,  and  3.14  are  derived  in  Appendix  D.  For  the  first 
order  finite  difference  equation. 


x    =    .i-f)±i*  Vj-J  £-j  a1 


W|  4   ;,-,  -,      i        y2     1  ,~2 


(     '+£*.) 


(4.5) 
where    CO   -  COA  t  • 

A  discussion  of  CO  ^s  &n  order  here.  The  quantity  CO  is 
the  natural  frequency  of  the  oscillator  in  radians  per  second. 
Therefore, CO /\f  or  CO  has  the  units  of  radians  per  time  interval  and 
is  a  measure  of  the  resolution  of  the  solution.  A  small  value  of 
60  would  mean  that  the  system  is  being  sampled  many  times  per  cycle. 
For  example,  a  value  of  CO  equal  to  0.5  would  mean  the  system  was 
being  sampled  every  0.5  radians,  or  over  12  times  per  cycle.   A  value  of 
CO  over  1.0  would  mean  that  the  solution  has  a  very  poor  resolution, 
and  the  truncation  error  would  be  high. 

Equation  4.1  can  be  written  in  the  form 


or 


y     =  A  e  a  PAf) 


Hi  =A[ 


eaAt 


J 

(4.6) 
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Comparing  Eqs.  4.4  and  4.6  reveals  that 

e aAf    vs        X 

is  a  proper  comparison.  However,  Eqs.  4.3  and  4.5  show  that  both 

&  and  \  are  complex  numbers  for  the  underdamped  system.  There- 
fore, 

aAf  -       OjAt    ±i  a_  At 
e  e  e  (4  ?) 

where 

ar  =  -cj£  ;    a-  =  aW/-^ 

and  for  the  first  order  difference  scheme 

X   =    X  ±  /  Z 

where  _ 

c  ~~     /  '  -l  y  -.  \  and      c 


The  solution   \  can  also  be  given  in  the  form 
where 


r^Sjx**  zz  and  Q    =APC   SIN   -f- 

A  further  simplification  can  be  made  by  replacing  r      with 
e'K.     Thus, 

k  =  -   In    r 

and  )^   becomes     w  +  \  f)  (4.8) 

=  e  e 

Examining  Eq.  4.7  and  4.8  reveals  that 

a  A  t        vs       -  k 

r 

and 

a.  At       vs        Q 

are  proper  comparisons.  Let  the  relative  error  in  both  K      and  (y 
be  defined  by 

Applying  the  value  of  Q.  •   and  Q.    to  LL  anc*  1/    »    leads  to 

U-   %  and    V  =     —^ (4.9) 
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The  relative  error  terms,  n  and  ]/     ,    provide  an  indication 
of  the  accuracy  in  the  numerical  solution.  This  can  be  seen  by  com- 
paring the  complete  solutions.   Substituting  Eq.  4.3  into  Eq.  4.1 
yields  the  exact  solution 

y  =  A   „[-£*»   ±  J  u^/pr^]  j /^t 

(4.10) 

since t=tJ/\f°      Substituting  Eq.  4.8  into  Eq.  4.4  yields  the  solu- 
tions for  the  finite  difference  approximations 

\-K    ±    i  6}  J 
y  =  A   eL 

(4.H) 

Applying  the  relative  error  terms  U,     and  ]/      given  by  Eq.  4.9  to 
Eq.  4.11,  the  finite  difference  solutions  take  the  form 

(4.12) 
Thus,  a  comparison  of  Eq.  4.12  with  4.10  shows  that  if  u      and"V 
are  both  equal  to  1.0,  Eq.  4.12  is  the  exact  solution. 

Examining  Eq.  4.12  reveals  that  IX     and  ]/      have  a  specific 
role  in  the  error  analysis.   In  Eq.  4.12,  —  Z/LCJis  a  measure  of  the 
rate  of  attenuation  in  the  finite  difference  solution.  Likewise, 
LLU)J  I  _}"Jis  a  measure  of  the  frequency  of  the  oscillations  in  the 
finite  difference  solution.  Therefore,  ]/     will  be  an  indication  of 
the  error  in  attenuation,  and  LL  will  be  an  indication  of  the  error 
in  frequency  of  the  numerical  sol 


— K 

~       For  the  overdamped  problem  C        does  not  change.   However, 

<?     becomes  an  attenuating  factor  and  Eq.  4.12  becomes 


[-2/£cu  ±  /^^\/^T/~JjAt 


y  =  A  e 


where  p    is  given  in  Eq.  4.9,  and  LL 

/  9 


M=  -?. 


CO 


V^T 


(4.13) 
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The  expressions  for  \  ,  ^    ,  Q    ,  and  /A  for  Eqs,  3,4S  3,10,  and 
3.14  are  derived  in  Appendix  D.  Table  I  lists  the  values  of  \    , 

K   ,  and  Q   ,  and  Table  II  lists  the  values  for  / Q    for  all  three 
numerical  approximations.   The  subscripts  c9  1,  and  2  refer  to  conven- 
tional, MOD  1,  and  MOD  2  respectively. 

In  conclusion,  a  proper  comparison  of  the  solutions  to  the 
three  finite  difference  approximations  with  the  exaet  solution  has 
been  shown  to  take  the  form 

_k ..  6 

]/   —  — r and 


for  the  underdamped  problem  and 

„_    K  _      _  =  g 


I 


for  the  overdamped  problem,  A  computer  program  thatPmakes  the 
necessary  comparisons  was  written  in  the  FORTRAN  IV  language  for  the 
IBM  360  digital  computer.  The  results  obtained  from  this  program 
follow. 


TABLE 
LIST  OF  THE 


«-«* 


Kc=r/nl   ''^° 


I    +  £Q 


.2 


fl I    =    4  PC    S  /  A/    *  I 


2  f  —2  I     —+  \2 

-  b( 


"2 


X,   = 


(/  -  -fu/osV/-  £*-£ 


u3a 


j-  /   .        >-  -  w 


H 


(/-   £*  +    T  ) 

K,    —  t  ln\      ; z-= 2]5— \ 

2_     Y2— I      j_      -4 
Q     =     ARC      SIN     1    I      ,    ~    ;>-:,,      Jj       "       -« 

i/+  c«*  Sr  ) 

*      Ail  ir  *'*■  7k  *'- £<*)* 


^         cu   -   Q  u)    -    T  co 


fc=4*C    S//V  J    [-^  ^Xfl^,^ 


6      ^     '    /?* 


25 


TABLE  II 
LIST  OF  THE  OVERDAMPED  SOLUTIONS 


_2y2  _Z  ,      -  *     \Y2L 

CO 


C 


'ft  =  ^c  S/NH  -  ,  _  »& 


i  Qt    —ARC     SINH 


CO     ^      -       CO    +/2"    CO 


par^i-  sr*j& 


'/i  co  C,    -   co     +  tt  co 
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CHAPTER  V 
RESULTS  AND  CONCLUSIONS 

Results 

Figure  2  is  a  plot  of  the  data  in  Table  III,   For  this  plot, 
r    is  equal  to  0.0,  and  the  system  is  undamped.  The  undamped  solu- 
tion will  have  no  attenuation  and  therefore,  the  V    term  is  meaning- 
less. 

With  Y     =0.0,  the  differential  equation  becomes 

m  y"+  k  y    =  q 

(5.1) 
Recall  that  the  MOD  2  method  was  derived  by  maintaining  the  Hermit ian 
form  of  the  second  derivative.   As  a  result ,  the  MOD  2  method  pro- 
duces very  little  error  for  all  values  of  CO  .  Since  the  first  de- 
rivative is  not  present  in  Eq«  5-1,  the  MOD  1  method  can  be  expected 
to  have  more  error  than  MOD  2.  However ,  it  is  found  that  the  error 
in  MOD  1  is  approximately  the  same  as  the  error  associated  with  the 
conventional  method. 

The  significant  gain  in  accuracy  by  the  employment  of  MOD  2 
means  a  larger  time  interval  could  be  used  to  approximate  the  system 
without  loss  of  accuracy.  For  exampl@s  if  the  natural  frequency  of 
the  system  equaled  180  radians  per  second,  a  time  interval  of  0.005 
seconds  would  describe  the  system  with  0,2  per  cent  error  in 
using  MOD  2.  To  achieve  this  accuracy  with  the  conventional  method, 
a  time  interval  of  0.0015  seconds  must  be  used.  This  means  that  the 
time  interval  can  be  increased  by  a  factor  of  three  by  using  the 
MOD  2  approximation  for  this  problem. 

As  the  damping  is  increased ,  the  coefficient  of  the  first 
derivative  increases s  and  the  first  derivative  begins  to  have  a  sig- 
nificant influence  upon  the  response.  The  MOD  1  method  should  become 
more  accurate  as  f      increases  since  it  is  derived  by  maintaining  the 
Hermitian  form  of  the  first  derivative.  Figures  3  and  4  are  plots  of 
the  data  presented  in  Table  I?  with  [    equal  to  0.9,  Since  £   is 
less  than  1.0,  the  underdamped  comparisons  are  used.   As  anticipated, 
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the  MOD  1  approximation  is  seen  to  be  the  most  accurate.   For  example, 
in  Figure  3  at  a  resolution  of  0.8  the  MOD  1  approximation  has  3  per 
cent  error  in  Q      .   For  this  same  resolution,  the  conventional  approxi- 
mation has  42  per  cent  relative  in  Q 

Figures  5  and  6  are  plots  of  the  relative  error  versus  f 
for  a  constant  resolution  of  0.5  as  listed  in  Table  V.   Figure  5 
shows  that  one  of  the  two  Hermitian  approximations  always  has  a  better 
value  of  LI  than  the  conventional  approximation.   For  low  values 
of  (_     the  MOD  2  approximation  is  best,  and  for  high  values  of  f     the 
MOD  1  approximation  is  best.   This  means  that  throughout  the  range 
of  r     ,  greater  accuracy  in  the  frequency  can  be  obtained  by  using  one 
of  the  Hermitian  approximations. 

The  error  in  attenuation,  ]/   ,  is  shown  in  Figure  6.  Note  that 
for  the  lightly  damped  system,  where  the  attenuation  is  small,  the 
conventional  approximation  has  the  least  error.  However,  this  does  not 
imply  that  the  total  solution  given  by  Eq.  4.12  is  most  accurate  when 
the  conventional  approximation  is  employed  since  Figure  5  shows  that 
the  conventional  approximation  has  the  largest  error  in  LL     for  the 
lightly  damped  system.  For  values  of   T  between  0.45  and  0.65,  the 
MOD  2  approximation  is  best,  and  for  values  of  f     greater  than  0.65, 
the  MOD  1  approximation  is  best.   Thus,  Figure  6  shows  that  the  Her- 
mitian approximations  have  less  error  in  attenuation  for  values  of  I 
greater  than  0.45. 

Analysis  was  not  made  for  the  critically  damped  oscillator. 
However,  as  T     approaches  1.0  from  either  side,  Figure  5  shows  that 
LL   becomes  very  large.   In  order  to  understand  what  happens  to  the 
numerical  solution  for  ,9  <C,  T^-Lls  it  is  necessary  to  make  other 
comparisons.   For  example,  the  exact  solution,  Eq.  4.10,  and  the 
solution  for  the  finite  difference  approximations,  Eq.  4.12,  can  both 
be  given  in  the  following  form 

y,A  ;[-£  i;  VT^jca^ 


and 

y--A 


f-O*  '>vp~rQ  U)AfJ 
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for  F  <C  1.0.  Thus  a  plot  ofv//  —  T2  &n&i£\jl—  yx  versus  T  shows 
the  error  in  the  oscillatory  part  of  the  solution  for  the  underdamped 
system.  For  the  overdamped  system,  f   ^  1,0,  Eqs.  4.10  and  4.12 

become       ["  £  ±  VF^  /  ]  U)Af  J 

Thus,  for  this  range  of  £   a  plot  of  /  V-L,  /  an£*  U.  v/ ^—  /  versus  L 
will  be  the  equivalent  comparison  for  the  overdamped  system. 

Examination  of  the  equations  in  Table  I  reveals  that  the 
numerical  solution  will  be  critically  damped  for  the  values  of  f slightly 
less  than  1.0.  Therefores  to  show  the  critical  damping  point  for  each 
finite  difference  approximation  it  is  necessary  to  use  the  comparisons 
for  the  overdamped  system  for  values  of .95^TjCl.O.   Equation  4.13 
defines  LI   for  the  overdamped  system  as 

u=  J4 

therefore,  b 

fO 


b  w  (5.2) 

Using  Eq.  5.2,  it  is  possible  to  evaluate  the  overdamped  comparisons 
for  values  of  C   <C  1.0. 

Figure  7  is  a  plot  o£n\Jj  -  Hor  UyJ'f2—   /for  values  of  L, 

between  0.9  and  1.1  with  Gl)  -  0.5.  The  MOB  1  approximation  has  its 
critical  damping  point  closest  to  1.0,  and  the  conventional  approxima- 
tion has  its  critical  damping  point  furthest  from  1.0.   As  U)     is  de- 
creased and  the  system  is  sampled  more  often,  the  critical  damping 
points  of  all  the  numerical  solutions  will  approach  1.0.   However,  the 
Hermitian  approximations  will  always  be  critically  damped  closer  to 
1.0  than  the  conventional  approximation.   The  fact  that1 each  numerical 
solution  has  its  own  critical  damping  point  that  is  less  than  1.0 
explains  the  large  values  of  IX     in  this  region. 

For  the  overdamped  system  a  problem  arises  when  f  r\  is  greater 
than  1.0.  Examination  of  the  equations  in  Tables  1  and  II  reveals  an 
unstable  zone  for  these  finite  difference  equations  at  this  point. 
Although  f    U)      -  1.0  is  essentially  the  limit  for  numerical  stability, 
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the  accuracy  of  the  difference  approximations  becomes  very  poor 
before  the  unstable  zone  is  reached.  As  the  stability  limit  is 
approached,  both   LL  and  1/     approach  infinity.  Figure  8  shows  the 
stability  limits  on  all  three  finite  difference  approximations.  This 
plot  was  made  using  the  equations  in  Tables  1  and  II. 

Conclusions 

From  the  results  presented  heres  the  conclusion  can  be  made 
that  the  Hermitian  equations  can  be  used  to  obtain  more  accurate 
numerical  solutions  for  the  same  size  increment  than  can  be  obtained 
by  first  order  central  finite  difference  equations  when  applied  to 
differential  equations  containing  derivatives  of  more  than  one  degree. 
Further,  this  Increase  in  accuracy  is  not  accompanied  either  by  an 
increase  in  complexity  of  the  resulting  matrix  when  applied  to  boundary 
value  problems  or  by  a  significant  increase  in  computation  time  when 
applied  to  initial  value  problems.  Consequently,  use  of  the  Hermitian 
equations  would  be  most  advantageous  when  computer  core  storage  space 
is  limited  and  accurate  results  are  required. 
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APPENDIX  A 

THE  DERIVATION  OF  THE  CONVENTIONAL  EQUATIONS 
FOR  THE  FIRST  AND  SECOND  DERIVATIVE 

The  first  order  central  finite  difference  equations  can  be  de- 
rived from  the  Taylor  series.  The  stations  are  lableled  station 
j+1,  station  j,  and  station  j-1.  The  interval  between  stations  j+1 
and  j  is  called  A 

The  Taylor  series  for  y..,  is 

//  a    '       tt      "       A^  "'     A    a 


r 
A   2" 

~y,    + 

si  V 


*  »   » 


(a.i: 


For  y.  .,  the  Taylor  series  becomes 


A^  //    A3   '" 


Hj-*  =  $-   A^  +  -yj-  -jm 


A*"  m  a5"  ff- 

?!    V         S!  V  (a. 2: 


Subtracting  Eq.  A.I  from  Eq.  A. 2  gives 


'         A3  "/   A6  z 


Solving  for  y'.  leads  to 


(A. 3) 


V      2& 


\        A-2  M        A*   2" 


This  is  the  first  order  central  finite  difference  equation  for  the 
first  derivative. 

Adding  Eqs.  A.I  and  A. 2  gives 

_  £!y*   +  .  -  • 

360^  (A. 5) 

Equation  «\.5  is  the  first  order  central  finite  difference  e< 

for  the  second  derivative. 
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The  truncation  error  in  each  equation  is  usually  approximated  by 
the  first  term  in  the  Taylor  series  that  is  not  used  in  the  difference 
approximation.  For  example,  the  error  in  Eq.  A.  4  is 


A2 


y, 


6     v 

and  the  error  in  Eq.  A. 5  is 

12.      *J 
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APPENDIX  B 

THE  DERIVATION  OF  THE  HERMITIAN  EQUATIONS 
FOR  THE  FIRST  AND  SECOND  DERIVATIVES 

Taking  the  derivative  of  Eq.  A.l  in  Appendix  A  gives 

Likewise  Eq.  A .2  becomes 

/  /  //    A2  ///     A3  i^ 


*   *    ' 


Adding  Eqs.  B.l  and  B.2  gives 

/        /         /      ±      '"  A*     * 

Dividing  Eq.  B.3  by  6  and  solving  for   A   y' ' '  leads  to 


6 

A^"77-   ±1.,'   _,„'.  u'    \-&uZ. 

6    J 


yj  =  m..-^*  Hhirfhi* 


From  Eq.  A. 4  of  Appendix  A, 

6  ^J      a£S  -y 


#  "Wl%"     ^"' 


-^;-  *;**, 


/  < » 


Thus,  substituting  Eq.  B.4  into  Eq.  A. 4  gives 

yJ=ikk/H-  yj-H&-iyj'*  yj'-<) 


(B.l) 


A  *  ^   A5"  J3 


(B.3) 


(B.4) 


(A. 4) 


(B.5) 
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Solving   for  y     .    -  y.   ,    gives 


u       -  u        -A 


yJ- 


1°     V  (B.6) 

Equation  B.6  is  the  Hermit ian  equation  for  the  first  derivative. 


The  second  derivative  of  equation  A.l  is 
//        //         "/   A2  &       A3. ,  % 

*jh   =yJ+  A  y'  f^y'  +    *y' 
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A*  .« 


y.  ^Ti^y, 


■     6  ^ 


Ar 


+.*.•• 


Likewise  Eq.  A.  2  becomes 

//  ,/  ft*      Az  i£     a3    -2" 


,J6 


A   Zl  A 

XfXJ       l^o^J 

Adding  Eqs.  B.7  and  B.8  gives 


yn     v 


Dividing  Eq.  B.9  by  12  and  solving  for  A   y    yields 

12 


M,^-_   / 


// 


v  I  //        //       //  \ 


A*,^ 


/z  -£/       /z  j-t/*/    ~f 

From  Eq.  A. 5  of  Appendix  A, 


/?*  -c; 


+■  •  *  * 


A^..^ 


IM     M  ^J 


^/    A2(#"  *%* 

J  GO  Jj 
Substituting  Eq.  B.10  into  Eq.  A. 5  gives 

yS*Myj«-zyj  *  yj-\-Ayj«-*yj"+  yj'-l 

+    —   w        +  •  *  * 


(B.7) 


(B.8) 


(B.9) 


(B.10) 


(A. 5) 


(B.ll) 
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Then,  solving  for  (y  .  -  2y.  y-  y.  ,)  gives 
I  \    AZ       "       ft 

\u     - 2y    +  y    \=-£-  u 


y.    +ioy    *  y 

A4  JZL 

—  y 


(B.12) 
Equation  B.12  is  the  Hermitian  equation  for  the  second  derivative. 
A  comparison  of  the  conventional  equations  in  Appendix  A  with 
the  Hermitian  equations  derived  here  reveals  that  the  Hermitian 
equations  have  less  truncation  error.  This  is  to  be  expected  since 
the  error  term  in  each  first  order  equation  was  used  to  derive  the 
respective  Hermitian  equation. 
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APPENDIX  C 

DISCUSSION  OF  THE  TRUNCATION  ERROR 
IN  THE  MODIFIED  HERMITIAN  EQUATION 

In  the  MOD  2  method  of  analysis,  the  Hermitian  equation  for 
the  first  derivative  is  modified.  Equation  B.6  is  the  Hermitian 
equation  for  the  first  derivative 


(B.6) 


This  equation  can  be  written 


sllkr^) -£&,*#*  */-/] 


A  £  % 
360  V  (c.i) 


The  conventional  equation  for  y' .  can  be  written 


\y^r  y 


'         A2 


H    *1tH/ 


/// 


4L  VJ+i       V-'. 

Addition  of  equations  CI  and  C.2  gives 


(C.2) 


C 


JZA 


y,  -  y,  \=  ir-\y! 


T     /I      ^J       360*J 


(C3) 

Equation  C.3  is  the  modified  Hermitian  equation  for  the  first 
derivative  that  was  used  in  the  MOD  2  approximation.  The  truncation 
error  in  this  equation  is 

A^        /" 

/x      Vj 

This  truncation  error  is  one-half  the  truncation  error  in  the  first 
order  equation  for  the  first  derivative.  Therefore -t    in  the  MOD  2 
approximation  to  the  differential  equation,  there  is  very  little  gain 
in  accuracy  associated  with  the  first  derivative. 
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Likewise,  for  the  MOD  1  approximation,  the  Hermitian  equation 
for  the  second  derivative  must  be  modified.  The  modified  equation  is 

-   J±-  u —  u 


(C4) 


The  truncation  error, 


A*  & 


y 


12    v 

is  equal  to  the  truncation  error  of  the  first  order  difference 
equation  for  the  second  derivative.  Therefore,  in  the  MOD  1  approxi- 
mation, there  is  no  gain  in  accuracy  associated  with  the  second 
derivative. 
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APPENDIX  D 

SOLUTIONS  FOR  THE 
FINITE  DIFFERENCE  APPROXIMATIONS 


The  conventional  approximation  will  be  used  as  an  example  for 
this  appendix.  The  conventional  finite  difference  expression  for 
the  differential  equation  is 


m 


&t\H)+r3-yJ  *yj-, 


+2At(  yJ+l       yJ~' 


+  ky,   =   0 


(D.l) 


When  the  assumed  solution 


yu  =  A\ 


j 


is  substituted  into  Eq„  D.l,  the  result  is 

m 


I&K-a-'CK-irK- £)♦.*]  **• 


.(At? 

=  o 

For  an  arbitrary  solution, 


mi  4  -^ )  *  c 


(A 


2ht 


Xc"X      ]+  kXc    =  ° 


or 


liAtf  "5AtJ 


m        c 


fc 


(A^J 


K* 


m 


1^*7     2AfJ 
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By  using  the  quadratic  equation 


X  = 


l.  -2HL 
k     (At 


jj  -     V  [k  ~  (^~  j£f+ 2af]  |^2.  J^J 


2  L  f At)z 


+ 


2AtJ 
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or 


K= 


Zm  -k(AA±~\J\zm  -k  (Atf|-V[m-  -jc<£tf    ] 


2 


m 


f&t 


When  the  substitutions 


U)  ~ 


-to- 
rn 


,  c- 


2tnoj 


and    co-cJAf 


are  made,   \   becomes 


i-icf 


isv^ 


/  *  ^ujZ 


I    +-£55 

To  make   X   complex,  take  ~\T^   outside  the  radical.   Then 


-i 


-2 

-*-  CO 


±/oj 


I  .-.*• 


K- 

The  solution    X   is  now  of  the  form 
'v  c 

Xc=.  x  +  i  z 

To  get   X   into  the  form 

the  following  relationships  are  used. 


/   -   C    "fw 


(D.2) 


5  =  Vl?773 

Therefore, 


and 


9c-arc  s/n-| 


_+ll;-*aJ]+-o'[/-S- 

<c    — ■— - — 


*>  ]  l* 


i  +  Cw 
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or 


6  =t-\/T^X^ 


'     +       £" 


«.  K    =  -  /n  r, 


then 


Similarly, 


Kc-  -  In  +  -^  //   -   /   CO 

V   I    +    [ft  (D.3) 


Y         /  -  £*coz 


(D.4) 
For   the  overdamped   system,    i     Q      must  be   found.     From  Eq.   D.4, 


sin  $e  -  t    /  sj  -  cv  - 


*-     CO 


*,7»2 


Then,  '  o 


_  *       AooV  £coV£co* 


V    /   —      c,  oo  c 

since 

/     S/A/  ^     -   SINH    I   Q 
Therefore, 


/  $     -ARC     SI  NH  t         /£gg    -    Of2   +^£):2' 

/  -  e^ 

(D.5) 
The  equations  for  \  ,  K"  ,  £  ,  and  i  £)  for  MOD  1  and  MOD  2 
are  derived  in  a  similar  manner. 
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